Dispersive bottleneck delaying thermalization of turbulent Bose-Einstein Condensates 
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A new mechanism of thermalization involving a direct energy cascade is obtained in the truncated 
Gross-Pitaevskii dynamics. A long transient with partial thermalization at small-scales is observed 
before the system reaches equilibrium. Vortices are found to disappear as a prelude to final thermal- 
ization. A bottleneck that produces spontaneous effective self-truncation and delays thermalization 
is characterized when large dispersive effects are present at the truncation wavenumber. Order of 
magnitude estimates indicate that self-truncation takes place in turbulent Bose-Einstein conden- 
sates. This effect should also be present in classical hydrodynamics and models of turbulence. 
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The Gross-Pitaevskii equation (GPE) furnishes a dy- 
namical description of superfluids and Bose-Einstein 
Condensates (BEC) that is valid at very low tempera- 
tures [Tj. The GPE dynamics is known to produce an 
energy cascade that leads to a Kolmogorov regime of tur- 
bulence [2j [3] . Such turbulent regimes were studied in 
low-temperature experiments in superfluid 4 He [4 and 
in BEC [5 . The truncated GPE (TGPE), obtained by 
performing a truncation of Fourier modes, can also de- 
scribe the (classical) thermodynamic equilibrium of ho- 
mogeneous BEC [6]. The TGPE (microcanonical) equi- 
librium is known to present (when varying the energy) a 
condensation transition [6j [7] . 

In the context of classical fluids, the (conservative) 
truncated Euler dynamics is known to possess long- 
lasting transients describing dissipative phenomena |5]. 
With this motivation, we study in this Letter the TGPE 
thermalization that arises from the GPE turbulent en- 
ergy cascade. Here is a short summary of our results. 
Partial thermalization is observed at small-scales during 
a long transient regime; vortex lines then disappear and 
final thermalization sets in. A bottleneck that delays the 
final thermalization is characterized when large disper- 
sive effects are present at truncation wavenumber. 

The TGPE is obtained from the GPE describing a ho- 
mogeneous BEC of volume V by truncating the Fourier 
transform of the wavefunction ip: ip^ = for |k| > /c max 
[TJ [6] . Introducing the Galerkin projector Vq that reads 
in Fourier space Vg[^\z] = 0(k max — |k|)^k with #(•) the 
Heaviside function, the TGPE explicitly reads 
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where m is the mass of the condensed particles 
and g = 47rah 2 /m, with a the s- wave scatter- 
ing length. Madelung's transformation ip(x,t) = 

• exp [i^0(x, £)] relates the (complex) wavefunc- 



p(x,t) 



tion ip to a superfluid of density p(x,t) and velocity 
v = V0, where h/m is the Onsager-Feynman quantum 
of velocity circulation around the ip = vortex lines pQ . 



When Eq.Q is linearized around a constant ip = ipo, 

- 1/2 

the sound velocity is given by c = (g\ipo\ 2 /m) with 
dispersive effects taking place at length scales smaller 

1/2 

than the coherence length £ = (h 2 /2m\ip \ 2 g) that 
also corresponds to the vortex core size. In the TGPE 
numerical simulations presented in this Letter the den- 
sity p = mN/V is fixed to 1 and the physical constants in 
eqs.ffl are determined by the values of £/c max and c = 2. 
The quantum of circulation h/m has the value c^/y/2 
andV = (27r) 3 . 

Equation ([!]) exactly conserves the energy H = 
Jd 3 x (^IWI 2 + § [PgM 2 ] 2 ) and the number of par- 
ticles N = jd 3 x\ip\ 2 . Using Fourier pseudo-spectral 
methods the momentum P = y J d 3 x (ip\7ip — ip\7ip) is 
also conserved with dealiasing performed by the 2/3-rule 
(tax = 2/3 x M/2 [9] at resolution M). 

We now study the thermalization of the superfluid 
Taylor- Green (TG) vortex, a flow which develops from a 
spatially-symmetric initial condition prepared by a min- 
imization procedure. The TGPE integrations are per- 
formed with a dedicated pseudo-spectral code that uses 
the symmetries to speed up computations (see reference 
[2]). Up to 512 3 collocation points are used and the co- 
herence length is determined by £/c max = 1.48. 

Vortices and density fluctuations are visualized on Fig. 
[I] The short time behavior, see Fig. [I]a-b, corresponds 
to the GPE superfluid turbulent regime previously stud- 
ied in [2 J. A new TGPE thermalization regime where 
vortices first reconnect into simpler structures and then 
decrease in size with the emergence of a thermal cloud is 
present at latter times, see Fig{T]c-d. 

To further study this relaxation dynamics, we express 
the energy per unit volume E tot = (H — mc 2 N)/V + m ^- 
as the sum of three (space-averaged) parts [2 j: the 
kinetic energy E^ n = (l/2(y/pVj) 2 }, the internal en- 
ergy Ei nt = ((c 2 /2)(p — l) 2 ) and the quantum energy 
E q = (c 2 £ 2 (dj s/p) 2 ) • Parseval's theorem allows to de- 
fine corresponding energy spectra: e.g. the kinetic en- 
ergy spectrum E^ n (k) as the (solid angle integral) of 
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FIG. 1: 3D visualization of density at t = 5, 10, 31 and 55 at 
resolution 512 3 (see below Fig|2Jb). Vortices are displayed as 
red isosurfaces and clouds correspond to density fluctuations. 
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E^in(k) can be further de- 
composed into compressible E^ in (k) and incompressible 

E Ln( k ) P arts > usin s (Vpvj) = (Vp v j) c + (Vp v jY with 

V • (s/pvj) 1 = 0. The temporal evolution of i^kin, ^kin> 
^kin' ^q+^int is displayed in FigJ2la and the correspond- 
ing energy spectra on Figj2]c-d. 

Three evolution phases are apparent on FigJ2]a. The 
first phase, for t < 15, corresponds to the GPE regime 
previously studied in [2 1. In the second phase, for 20 < 
t < 70, the spectral convergence of the GP partial dif- 
ferential equation is lost and the dynamics is influenced 
by /c max . Partial thermalization starts to take place at 
large wavenumbers where £^kin(^) ~ k 2 (see FigJ2]c). Fig- 
ure [2]b shows that this phase is delayed when the res- 
olution is increased at constant £/c max . When t > 80 
the system reaches the thermodynamic equilibrium with 
equipartition of energy between E^ in and E q + E[ nt , see 
Fig. J2Jd. Finally, E^ in vanishes before final thermaliza- 
tion (see FigJ2a-b); the total disappearance of vortices is 
observed on corresponding 3D visualizations (data not 
shown). Similar relaxation mechanisms are also present 
in models of hydrodynamic turbulence and the truncated 
Euler dynamics [8] [TDJ [11] . 

We now focus on a characterization of thermodynamic 
equilibrium that will allow us to account for the absence 
of vortices and the equipartition of energy in the final 




FIG. 2: a) Temporal evolution of energies E^ in , E^ in , E^[ n 
and .Eq + E[ nt . At large times, the incompressible energy 
vanishes and equipartition of energy between E^[ n and Eq + 
E[ nt is observed. Resolution of 256 3 . b) Temporal evolution of 
El in at resolution of 64 3 , 128 3 , 256 3 and 512 3 with constant 
^max = 1.48. c-d) Energy spectra at t = 55 and t — 11 
resolution 512 3 and 256 3 respectivelly. 



states. Microcanonical equilibrium states are known to 
result from long-time integration of TGPE [6] [7] . Grand 
canonical equilibrium states are given by the probabil- 
ity distribution P st [^] = Z~ x exp[-/3(H - fiN)}. They 
allow to directly control the temperature (instead of the 
energy in a microcanonical framework). These states can 
be efficiently obtained by constructing a stochastic pro- 
cess that converges to a realization with the probability 
PstM [12] • This process is defined by a Langevin equa- 
tion consisting in a stochastic Ginbzurg-Landau equation 
(SGLE): 
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(2) 



where the white noise C( x ^) satisfies (C(x,£)£*(x', t')) = 
S(t — t f )5(-x — x'), j3 is the inverse temperature and /i the 
chemical potential. Using this algorithm in [T2] the mi- 
crocanonical and grand canonical ensembles were shown 
to be equivalent and the condensation transition reported 
in [6j [7] identified with the standard second order A- 
transition (see insets on Figj3J). Note that c p would be 
very difficult to obtain from microcanonical results. 

At low-temperature the partition function Z can be 
exactly computed by the steepest-descent method [T2] . 
In particular, the mean value of the condensate ampli- 

2 - . E AT f [_4m/x_ 

~ 9 V/3/i^0 



tude reads |^o 



where M = 



^Lx^/^ 2 i s the total number of modes and fo[z] = 
3z - 9z 3 / 2 cot" 1 (y/z) /4. 

In numerical simulations of Eq.([2|, /i is adjusted to fix 
the density (or the pressure p) [12]. The inverse tem- 
perature is normalized as /3 = JV/VT. With this choice 
of parametrization the A-transition temperature T\ is in- 
dependent of M. Data from SGLE and low-temperature 
calculation are confronted on Figj3]and seen to be in good 
agreement. 
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FIG. 3: Temperature dependence of the energies E^ in (stars), 
El in (diamonds), E^[ n (circles) and E^ + E^t (pentagrams) at 
constant density. Insets: (right) Temperature dependence of 
the condensate fraction \i/jo\ 2 /p; (left) Specific heat c p — -^ | 

at resolution 128 3 . 



The temperature dependence of the different energies 
is displayed on FigJ3l Observe that E^ in vanishes at tem- 
peratures T/T\ < 1/2. This explains the disappearance 
of vortices in Fig J2] above as the corresponding final tem- 
perature is well below T\ (see corresponding values of 
energies on Figj2ja). At low temperature equipartition 
of energy between E^ n and E q + £j nt is also apparent on 
Figj3] Note that a larger fc max implies, by equipartition, 
a lower temperature. The corresponding dissipative ef- 
fects are thus smaller explaining the thermalization delay 
apparent on FigJ2jb. 

We now turn to the study of dispersive effects on the 
thermalization of the TGPE dynamics. In order to in- 
vestigate dispersive effects, the TG initial condition de- 
scribed above is prepared at fixed £ = a/2/20 and varying 
resolution: 64 3 , 256 3 and 256 3 . The three initial condi- 
tion thus represent the same field at different resolutions. 

The evolutions of the energies of the three runs are 
shown on Fig(4ja. They are identical until t ~ 5 where 
the run at resolution 64 3 starts to lose spectral conver- 
gence. All runs appear to have completely thermalized at 
t « 20. However the kinetic energy spectra corresponding 
to this time, displayed on FigJ4]b, shows clear differences 
between runs. The high- wave number modes are thermal- 
ized in the 64 3 run but they decay at large-/c at higher 
resolutions. In the 256 3 run, two zones are clearly dis- 
tinguished: an intermediate thermalized range (with ap- 




FIG. 4: a) Evolution of energies at £ = a/2/20 and resolution 
64 3 , 128 3 and 256 3 . b) Energy spectrum E kin (k) at t = 19.8 
for the three TG runs, c) Evolution of E^[ n (k), solid red lines 
are fits of the form Ak 2 exp [— j 2 k 2 ] (see text), d) Evolution 
of k c . Curves i-iv: £ = 2>/2/5, K n = 4, E [n = 0.1, 0.2 , 0.4 , 1; 
v: £ = v^/lO, K n = 8, E in = 0.2; vi-viii: £ = v^/5, 
E in = 0.1, 0.2,0.4 (i-viii in resolution 64 3 ); ix-xi: Taylor- 
Green resolutions 64 3 , 128 3 and 256 3 . e) Parametric repre- 
sentation dk c /dt v.s. /c c //c max (same labels as d). 



proximative k 2 scaling) followed, well before /c max = 85, 
by a steep decay zone. 

The temporal evolution of E^[ n (k) for the 256 3 run dis- 
played in Fig(4jc. is well represented by a fit of the form 
A(t)k 2 exp[-j 2 (t)k 2 }, where A(t) and k c (t) ~ 7 _1 (t) <C 
fc max are increasing functions of t. Such a behavior of 
the energy spectra ensures spectral convergence and the 
dynamics is thus not influenced by /c max . Furthermore, 
we checked that A(t)k 3 (t) « E kin for t > 20 (data not 
shown). This new regime can be described as a (quasi) 
thermalization, with self-truncation at wavenumber k c 
and temperature T ~ E/k 3 , that spontaneously estab- 
lishes itself within the GP partial differential equation 
dynamics when the direct energy cascade is inhibited by 
a dispersive bottleneck for the energy transfer [17] . 

An open question is whether thermalization is simply 
delayed or, as in the Fermi- Past a- Ulam-Tsingu problem 
[13], completely inhibited in the self-truncation regime 
£fcmax —> oo- It is not feasible now to directly study this 
limit, within the TG framework, as it requires long runs 



at arbitrarily high resolution. To skip the TG cascade 
regime and directly study the self-truncated thermaliza- 
tion regime we use initial data generated by the SGLE 
instead of the TG vortex. To wit, we use Eq.(p]) with a 
variable truncation wavenumber /cj, n , set to a target value 
of fc c , smaller than the maximum truncation wavenum- 
ber /c max allowed by the resolution. This SGLE-generated 
initial data is then used to run the TGPE at a given value 
of ^k c with arbitrarily large values of £/c max . 

A number of such runs were performed at reso- 
lution 64 3 with various values of /cj, n , £, and ini- 
tial energy E in . The result of these computations 
are compared with the TG runs and displayed on 
Figj4]d. The self-truncation wavenumber is explicitly 
determined from E^ n by the integral formula k 2 = 
5 j^max ^E^ {n {k)dk/ J fcmax E kin (k)dk. A general growth 
in time of k c is apparent on Figj4jd for both the Taylor- 
Green runs and the SGLE-generated initial data, showing 
similar behavior. 

In order to check for a self-similar regime a paramet- 
ric Log-Log representation of dk c /dt v.s. k c is used on 
Figj4]e. With this representation, a self-similar evolu- 
tion k c (t) ~ t a corresponds to a line of slope (a — l)/a. 
Figure [4je shows transient self-similar evolutions termi- 
nated by a vertical asymptote corresponding to logarith- 
mic growth (a = 0). This self-truncation takes place at 
small values of k c /k max strongly suggesting that the self- 
truncation happens in a regime independent of cut-off. 
This regime should, in principle, be amenable to a de- 
scription in terms of wave turbulence theory along the 
lines of reference [15] . 

Let us concentrate now on estimations of order of 
magnitude relevant to physical weakly-interacting BEC. 
At very low-temperature, the GPE gives an accurate 
description of the (classical) dynamics of BEC pQ at 
scales larger than the mean inter-atomic particle dis- 
tance £ ~ |^o|" 2/3 , satisfying a <C £ <C £. At fi- 
nite temperature [TJ |6], Bogoluibov's dispersion rela- 
tion cj|(fc) = k 2 g\il) \ 2 /m + k 4 h 2 /4m and the rela- 
tion hujB(k eq ) = /3 _1 = k^T imply that phonons are 
in equipartition only for wavenumbers k < fc eq . The 
equipartition wavenumber thus satisfies £/c eq ~ jyj 1 * 
for T <C T* and £/c eq ~ ^T/T* for T > T*, with 
T* = £ 2 T x ji 2 and T A - h 2 /k B m£ 2 the condensation 
temperature (see ref. [1 ). Thus physical BEC at low- 
temperature have a natural cut-off for the equipartition 
range given by /c max = k eq (T). 

In experimental BEC the value of £/c eq is large because 
T*/T x ~ £ 2 H 2 < 1 and the corresponding TGPE has 
a large £/c max . This strongly suggests that, unless over- 
whelmed by other (non-TGPE) relaxation mechanisms, 
the thermalization slowdown caused by the dispersive 
bottleneck should be observable. It would thus be in- 
teresting to investigate in the future the suppression of 



aspect-ratio inversion observed in recent turbulent BEC 
experiments [5] in the context of the TGPE in a non- 
homogenous trapping geometry. 

In the context of classical hydrodynamics, self- 
truncation can take place in fluids with dispersion e.g. 
MHD flows with Alfven waves [12]. Models of turbu- 
lence of the Euler-a [16 type, where the r.h.s. of the 
incompressible Euler equation is multiplied by the opera- 
tor (1 — aV 2 ) -1 that penalizes the energy at small-scales, 
are also plausible candidates for self-truncation. 

In summary we investigated the thermalization dy- 
namics in turbulent BEC using the TGPE. Our main 
result is that a bottleneck delays the final thermalization 
when large dispersive effects are present at truncation 
wavenumber and produces an effective self-truncation. 
These effects are in principle observable in physical BEC. 

We acknowledge useful scientific discussions with G. 
During and S. Rica. The computations were carried out 
at IDRIS (CNRS) and visualizations used VAPOR [18]. 
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